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ABSTRACT 


The  problem  of  selecting  signal  processing  parameters,  particularly  detection 
thresholds,  so  as  to  optimize  downstream  tracking  performance  is  examined 
further.  Numerical  simulations  are  used  to  establish  the  validity  of  certain 
approximations  made  previously.  These  simulations  suggest  that  steady-state 
analysis  is  inadequate,  and  an  adaptive  threshold  optimization  scheme  is 
proposed  as  an  alternative.  Finally,  the  original  derivation  of  the 
Probabilistic  Data  Association  Filter  (PDAF) ,  upon  which  the  present  work  is 
founded,  is  augmented  to  account  for  finite  gate  size. 
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1.  INEROODCnOR 

A  critical  but  well-understood  issue  in  tracking  problems  involves  the 
inaccuracy  of  the  measurement  data,  which  is  typically  modeled  as  additive 
randan  noise  with  known  mean  and  covariance  [1-5] .  In  many  tracking  problems, 
particularly  those  arising  in  surveillance,  there  is  an  equally  critical  but 
less-understood  uncertainty  in  the  origin  of  the  received  data,  which  may  (or 
may  not)  include  measurements  from  the  target (s)  of  interest,  interfering 
targets,  or  random  clutter  (false  alarms) .  This  leads  to  the  problem  of  data 
association  or  data  correlation,  which  has  been  attacked  on  a  number  of  fronts 
[6-141  and  surveyed  in  [15-17] .  In  this  situation,  tracking  performance 
depends  not  only  upon  the  noise  covariances,  but  upon  the  amount  of 
uncertainty  in  measurement  origin.  In  some  of  the  approaches  cited  above 
[6-10] ,  this  dependence  is  explicit  and  is  characterized  in  terms  of  the 
detection  probability  PD  and  false  alarm  probability  PF. 

As  shown  in  Figure  1,  measurement  data  are  typically  provided  to  the  tracks  by 
some  sort  of  signal  processing  and  detection  algorithm,  where  the 
probabilities  of  detection  and  false  alarm  are  controlled  by  the  selection  of 
a  detection  threshold.  In  a  previous  stage  of  this  project  [20] ,  we 
established  a  quantitative  relationship  between  this  threshold  and  the  state 
error  covariance  matrix  in  the  tracker  downstream.  More  specifically,  it  was 
shown  that  if  the  tracking  is  done  with  an  extended  Kalman-Bucy  filter 
modified  to  use  probabilistic  data  association  (PDA)  [6-8,15] ,  then  its 
conditional  covariance  update  equation  (stochastic  Riccati  equation) ^ 

m 

*kik  -  ik ik-i  -  *  *<  (i) 

j*l 

can  be  approximated  by  the  deterministic  equation 


^The  notation  is  all  defined  in  [20] 
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Figure  1.  Tracking  system  block  diagram 
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*k|k  *  ^klk-l  “  <32(sk?pD'pF)Hks^iik 


(2) 


where  the  scalar  quantity  q o,  which  lies  between  0  and  1,  depends  upon  (via 
the  gate  volume  V  =  Cj^lfil  V2) '  *D'  and  PF.  Although  q2(Sj{;  pd'pf>  is  defined 
by  an  infinite  sum  involving  nested  integrals,  it  has  been  evaluated  using 
Monte  Carlo  integration  and  reduced  to  a  table  look-up  procedure  that 
facilitates  numerical  evaluation  of  (2) . 

In  [20],  (2)  was  iterated  to  a  steady-state  value  £(PD,PF)  for  a  particular 
tracking  example,  and  the  steady-state  RMS  position  error  I(PD,PF)  was 
displayed  on  a  contour  plot  with  coordinate  axes  PD  and  Pp,  called  a  tracker 
operating  characteristic  (TOC) .  A  typical  set  of  TOC  contours  is  shewn  in 
Figure  2,  and  a  corresponding  set  of  receiver  operating  characteristic  (ROC) 
curves  appears  in  Figure  3. 

Each  ROC  curve  represents  a  locus  of  possible  operating  points  for  the 
detector/receiver  that  is  providing  measurements  to  the  tracker:  a  particular 
setting  of  the  detection  threshold  selects  a  point  along  the  curve  and 
determines  the  values  of  PD  and  Pp  that  will  affect  the  tracking  performance 
via  (2).  Thus,  if  the  appropriate  ROC  curve  is  superimposed  on  the  TOC 
contours  [see  the  dashed  line  in  Figure  2] ,  the  dependence  of  tracking 
performance  on  detection  threshold  may  be  determined  graphically.  The  point 
marked  fl,  for  example,  represents  an  approximately  "optimal"  choice. 

The  purpose  of  this  paper  is  threefold.  First,  we  validate  the  procedures 
described  in  [20]  with  seme  more  extensive  numerical  simulations.  Second,  we 
propose  some  adaptive  detection  threshold  optimization  schemes.  Third,  we 
revise  the  PDAF  derivation  in  [€]  to  account  accurately  for  finite  gate  size. 

In  Section  2,  we  validate  the  procedures  described  in  [20]  by  comparing  the 
deterministic  approximation  to  the  stochastic  Riccati  equation,  the  PDAF- 
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Figure  2.  Tracker  operating  characteristic  (TOC)  contours 
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calculated  error  covariance  from  the  stochastic  Riccati  equation,  and  the  true 
PDAF  error  covariance.  Simulation  results  show  that  all  quantities  are 
comparable  for  a  reasonable  range  of  PD'PF*  Cutside  this  range,  PDAF 
performance  degrades  substantially  and  the  computed  quantities  are 
inconsistent,  especially  when  velocity  errors  are  considered.  Furthermore,  a 
filter  with  gains  based  on  the  deterministic  approximation  is  comparable  to 
the  PDAF  within  a  reasonable  range  of  PD'PF*  Finally,  we  give  several  reasons 
why  steady-state  values  of  the  deterministic  approximation  may  not  be  useful 
in  detection  threshold  optimization. 

In  Section  3,  we  propose  s<xne  adaptive  threshold  optimization  schemes  which  do 
not  require  iteration  of  the  deterministic  approximation  to  convergence.  In 
particular,  we  give  prior  and  posterior  algorithms  which  minimize  the  mean 
square  estimation  error  over  detection  thresholds  which  depend  on  observations 
up  to  the  previous  and  current  iteration,  respectively. 

In  Section  4,  we  revise  the  PDAF  derivation  from  [15]  to  account  for  finite 
gate  size.  It  turns  out  that  the  conditional  mean  update  is  the  same  as  the 
PDAF,  but  the  conditional  covariance  update  is  increased  by  an  additional 
term.  The  derivations  in  this  section  are  important  for  theoretical 
completeness  and  nay  be  of  practical  significance  under  certain  operating 
conditions. 
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2.  VALUATION  OP  APPROXIMATIONS 


Ihere  are  severed  reasons  for  obtaining  a  measure  of  FDAF  performance  as  a 
function  of  parameters  (e.g.,  PD,PF)  for  joint  detection/tracking  problems, 
including: 

1.  efficient  allocation  of  cormunication  and  computational  resources 

2.  optimization  of  parameters  for  inproved  performance 

A  computationally  efficient  method  of  obtaining  such  a  measure  has  been 
preposed  in  [20]  and  involves  iterating  a  deterministic  approximation  to  the 
stochastic  Riccati  equation.  We  shall  now  proceed  to  validate  this  procedure. 

Using  the  notation  of  [20],  we  let  be  the  output  of  the  FDAF,  f^ik  the 
PDAF-calculated  conditional  error  covariance,  £^|k  the  true  FDAF  conditional 
error  covariance,  and  the  output  of  the  deterministic  approximation  to 
the  stochastic  Riccati  equation,  i.e. , 


*k|k  *  *k|k-l  +  Hk2 

► 

Ik-1  =  E2k|k-llk-l  •  £0|0  5  E(*0}  = 

* 

m 

*kik  -  £kik-i  -  a-^0)HkSk^  + 

Sr  -flEkik-ifl'  +B  ;  Hk^lk-lSk1 

Ik— 1  *  ^^k-llk-l^'  +  222'  •  ^0|0  *  “  Eg 


(3) 


(4) 


•  y  —  T  »'*'■  'tr  .•*  •  •4* 
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*k|k  4  «»k|k*k|klYkl 

► 

|k-l  5  E{ik|k-lik|k-l|Yk”1} 

4|k  -  4 |k-l  -  '32^k'-pD-pFlKk44' 

4  4  fl^ik-ifl'  ♦  B  ;  s^44lk-l(Sk)'1 

*Sc|k-l  =  ^fc-llk-l*  *  ®G'  >  *010  *  *0 

Consider  the  standard  PDAF  assumption: 

PC^IY^-1)  -  JEk I k_!) 


(5) 


(7) 


where  N  denotes  a  normal  (Gaussian)  density.  If  this  assumption  is  satisfied, 
then 


^k |k  =  an^  E|<|k  *  ^k|k 


(8) 


To  validate  the  deterministic  approximation  (6)  to  the  stochastic  Riccati 
equation  (3),  we  compare  it  with  both  the  true  PDAF  error  covariance  and  the 
PDAF-calculated  error  covariance.  Since  the  standard  PDAF  assumption  is 
generally  not  satisfied,  the  PDAF-calculated  error  covariance  will  generally 
not  be  equal  to  the  true  PDAF  error  covariance.  In  the  above  notation,  we 
will  compare  £^|k  with  both  E{p£|k}  and  EtP^jk}* 
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2.1  Results 


All  results  are  for  the  numerical  problem  specified  in  [20] .  Monte  Carlo 
simulations  were  carried  out  to  compute  I^|k,  E{p£|k},  and  E{£^ |^}  for  the 
values  of  PD'PF  shown  in  Figure  4.  Sample  means  were  computed  from  ten 
(independent)  trials  using  random  variates  generated  by  IMSL  routines.  Hie 
values  M*E{pjj|k],  and  M“E{Plc j were  computed  for 


M  =  RMS  error  =  ,/trace  (P) 

M  =  RMS  position  error  Pn+P22 
M  =  error  volume  °C  V determinant  (£) 

M  =  position  error  volume "C  y determinant/  p^  pj_2\ 

\P21  P22/ 

In  Figures  5-12,  we  shew  plots  of  RMS  error  and  position  error  vs.  time.  From 
these  and  rhe  rest  of  the  data  we  observed  the  following  behavior: 

1.  For  Pd2  .4  and  Pp  ±  .06,  M*p£|k  =  M’EtP^^}  =  M*E{p£(k},  and 
typically  IM-Et^  j  , k  |  <  |M*  E{£^  j  k}-M'£^  |  k  | . 


2.  For  (PD,PF)  =  (.3,. 07)  and  those  metrics  M  which  reflect  velocity 
errors  (the  first  and  third) ,  PDAF  performance  degrades 

substantially,  and  M’£k|k'  M’E^|k^f  311(3  M>E^k|k*  are 

inconsistent. 

Because  the  deterministic  approximation  to  the  stochastic  Riccati  equation  is 
typically  comparable  to  the  PDAF-calculated  error  covariance,  a  question 
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2.2  Deterministic  Gains 


Before  we  present  simulation  results  for  the  PEAF  with  gains  based  on  P^|k-1' 
we  show  that  such  a  filter  is  optimal  (in  a  Bayesian  sense)  under  an 
assumption  stronger  than  the  standard  PEAF  assumption.  Let  ^  be  the  output 
of  the  PEAF  with  gains  based  on  |  k_1 ,  and  the  true  error  covariance  of 
this  filter,  i.e., 


*k|k  *  *k|k-l  + 

^k|k-l  =  ^k-1 1 k-1  ;  *0|0  =  a 


■2^tk  =  E^  (*k“*ktk)  ^"^klk^ '  ^ 

^kfk-1  5  E{ (Xk'Xkik-15  ^“^klk-15 '* 

Consider  the  following  assumption: 

Xj  and  Y^-1  are  jointly  Gaussian  for  j=l, . . . ,k. 
If  this  assumption  is  satisfied,  we  will  show  that 

pu*  I  l!*-1)  -  Httgnrf,  s£|k-i> 
and  consequently 


(10) 


(11) 
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*ic|k  *  El»k  I  **)  “4  £k|k  -  £k|k  112) 

Since  ij  and  Y^-1  jointly  Gaussian  for  j=l,...,k  implies  that  p(ij  |Y^-1)  is 
Gaussian  for  j=l, . . . ,k,  and  consequently 

P(£j  I  Yj-1)  ~  £j  |  j— 1> '  j=lf.../k,  (13) 

it  is  sufficient  to  show  that 

*k|k-l=4|k-l  ^  ^kjk-1  =  ^klk-l  (14) 

We  proceed  by  induction.  First  note  that  pf|0  =  £i|0*  Next,  for  any 
j=l , . . . ,  k-1  assume  that  £j|j-i  =  £j  j  j— x*  11111:115 

*jlj  -  4j-l  - 

=  Ejlj-i  -  q2(Sj;PDFPF)HjSjSi3  =  E{£j ,  j  |  Y^-1}  =  E{£jl;j}  (15) 

and  so 

Ej+llj  ***jlj£'  +fl“a'  =B(£j|j)E'  +SQS'  =  E(Ej+1|j)  “Ejrtlj  (») 

where  the  last  equality  follows  from  £j+1,  Y^  jointly  Gaussian  and  a  standard 
Kalman  filtering  argument.  Thus,  by  induction,  we  have  Ej|j_i  “  Pj  |  j— 1  for 
j-1, . . . ,k,  and  consequently  *  ^i^-i  and  *  ^k|k-l  38  recIuired* 
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2.3  Results 


In  Figures  13-20,  we  show  plots  of  M’2^  and  M*E{p£|k}  vs.  time  for  M  =  RMS 
error  and  position  error.  From  these  and  the  rest  of  the  data,  we  observed 
the  following  behavior: 

1.  For  PD  2  .5  and  Pp  1  .05,  =  ME{£j|k}. 

2.  For  (PD,Pp)  =  (.4, .06)  and  those  metrics  M  which  reflect  velocity 

errors  (the  first  and  third) ,  the  performance  of  the  PDAF  with  gains 
based  on  degrades  substantially,  but  normal  PDAF  performance 

does  not;  for  (PD,Pp)  -  (.3,. 07)  and  the  same  metrics,  the 

performance  of  both  filters  degrades  substantially. 


2.4  Steady  State 

In  the  ROC-TOC  approach  developed  in  [20] ,  the  deterministic  approximation  to 
the  stochastic  Riccati  equation  is  iterated  to  convergence.  We  oenote  this 
steady-state  value  as  £(PD,Pp) .  In  previous  simulations  [20] ,  we  have  shown 
that  £  exists  except  for  a  region  in  the  PD'PF  plane  where  the  deterministic 
approximation  to  the  stochastic  Riccati  equation  is  unstable  and  |k 
diverges.  However,  even  when  £  exists,  convergence  can  be  slew.  For  example, 
when  (PD,Pp)  *  (.6,.  04),  approximately  500  iterations  are  needed  for 
convergence  (see  Figure  21) .  Also,  both  existence  and  convergence  rate  are 
numerically  sensitive  to  the  initial  covariance  Pg|0.  For  these  reasons,  it 
is  questionable  whether  £  should  be  used  to  optimize  such  parameters  as  Pj^Pp, 
and  alternative  adaptive  approaches  are  developed  in  section  3. 
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2.5  Conclusions 

To  validate  the  deterministic  approximation  to  the  stochastic  Riccati 
equation,  we  compared  the  deterministic  approximation  with  both  the  true  PDAF 
error  covariance  and  the  PDAF-calculated  error  covariance.  Simulation  results 
show  that  all  quantities  are  comparable  for  a  reasonable  range  of  PD,Pp. 
Outside  this  range,  PDAF  performance1  degrades  significantly  and  the  computed 
quantities  are  inconsistent,  especially  when  velocity  errors  are  considered. 
Since  the  deterministic  approximation  is  typically  comparable  to  the  PDAF- 
calculated  error  covariance,  the  question  arises  as  to  the  performance  of  a 
PDAF  with  gains  based  on  the  deterministic  approximation.  We  first  showed 
that  such  a  filter  is  optimal  (in  a  Bayesian  sense)  under  an  assumption 
stronger  than  that  under  which  the  PDAF  is  derived.  To  evaluate  the  filter, 
we  compared  the  true  error  covariance  of  the  PDAF  with  gains  based  on  the 
deterministic  approximation  with  the  true  PDAF  error  covariance.  Simulation 
results  shew  that  these  quantities  are  comparable  for  a  reasonable  range  of 
PD'V  Finally,  we  considered  steady-state  issues.  In  view  of  the  slew 
convergence  of  £^|k,  it  appears  that  steady-state  values  of  the  deterministic 
approximation  should  not  be  used  to  optimize  such  parameters  as  PD'PF* 
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Figure  15.  H*Eg|k  and  M*E{E£|k}  for  (P^Pp)  *  (.8, .02) 
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3.  DETECTION  THRESHOLD  OPTIMIZATICN 

In  the  ROC-TOC  approach  developed  in  [20] ,  the  deterministic  approximation  to 
the  stochastic  Riccati  equation  is  iterated  to  convergence.  From  Section  2, 
we  denote  this  steady-state  value  as  £(PD,PF).  Then  some  metric  on  £  is 
optimized  over  the  values  of  PD'PF  which  satisfy  a  ROC  constraint.  There  are 
two  significant  problems  with  this  procedure  (see  Section  2) : 

1.  £(PD,Pp)  does  not  exist  for  certain  values  of  PD'PF 

2.  Even  when  £(PD,PF)  exists,  convergence  can  be  slow 

Thus  we  cure  led  to  examine  adaptive  approaches  to  the  problem  of  detection 
threshold  optimization.  The  adaptive  approaches  we  shall  consider  are  all 
time-dependent,  but  vary  in  data  dependence. 

Let  be  an  estimate  of  given  Y^-1,  and  £^[^-1  the  corresponding 

conditional  error  covariance.  We  make  the  standard  PDAF  assumption: 


P(Xr  I  Y^  ^)  -  N(*k|k-i'  £^  | ic— i ) 


(17) 


where  N  denotes  a  normal  (Gaussian)  density.  Thus,  if 


*k|k  =  2k|k-l  + 

m 

4 |k  -  *k|k-l  -  <l-*0>Hk£kH|!  +  Hk(Z*j2j2j-S2')H|i 

then 


(18) 
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*k|k  -  Ef%  I  ykl  «“>  £k|k  “  Ef*k|k*k|k  1  yk>  (M) 

Note  that  |k  depends  on  PD'PF  and  |  .  Since  is  a  function  of  Yk, 
PD,PF  can  (in  gene  rail)  be  taken  as  functions  of  Y]  for  j^k  but  not  for  j>k. 
In  the  sequel,  we  give  optimal  schemes  for  choosing  PD,Pp  as  functions  of  Y^-1 
and  Yk,  which  we  refer  to  as  prior  and  posterior  threshold  optimization, 
respectively. 

In  the  following  derivations,  we  ignore  certain  technical  questions  and  assume 
that  all  extrema  exist.  Let  Y^  be  the  space  of  realizations  of  Y^,  j=l,2, . . . , 
R  the  set  of  (PD,Pp)  which  satisfy  the  ROC  constraint,  and  (PD,Pp)  a  mapping 
from  Y^  into  R. 

3.1  Prior  Threshold  Optimization 

We  want  to  solve  the  following  functional  minimization  problem: 

(Pi)  minimize  E{  | k  | ^} 

over  (PD,PF) :  Yk_1  — »  R 

First,  we  have 

Etlj^-4^12}  *  E{tr  Ete^lY^1}}  (20) 

Since  E{^t||C|Yk”1}^0  V  Yk_1,  tr  E{£jt|k|Yk"1}i0  V  Yk_1,  and  consequently  (PI) 
reduces  to  the  pointwise  minimization: 
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(P2)  for  any  given  Y*4-1 

minimize  tr  E{£}<  uJY*-1} 
over  (PD,Pp)  €  R 

Let  (pQ/Pp)  be  the  optimal  value  of  (PD,Pp)  at  Y*-1  from  (P2) .  Then 

(P^Pp)  (Yk_1)  =  (P£,P;>  (21) 

is  optimal  for  (Pi) .  Now 

tr  Et^^lY14-1}  =  tr  -  q2(Sk?PDfPF)tr(Hk£ki^)  (22) 

Since  2.  Q.r  tr  (H^S^SJ^)  2.  0  r  and  since  q2  is  the  only  term  that  depends 

I.  i  “ 

on  PD,Pp,  tr  EfP^j^jY  }  will  be  minimized  when  q2 (S^/P^/Pp)  is  maximized. 
Thus  our  problem  becomes 

(P3)  for  any  given  Y^-^, 
maximize  q2 (£fc;PD,Pp) 
over  (PD,Pp)  €  R 

We  now  propose  an  algorithm  for  solving  (P3) .  First,  we  parameterize  PD,PF 
such  that 

R  =  {(pD(X),Pp(^))s  x  e  A}  (23) 
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where  A  is  a  closed  interval.  Since  g  corresponds  to  a  ROC  constraint,  we 
can  always  choose  x  to  be  PD,  PF,  or  the  detection  threshold.  Such  a 
parameterization  reduces  (P3)  to  a  line  minimization: 


(P4)  maximize  q2*(PD,PF) (x) 
over  X  6  A 


where  the  composition  *  is  defined  by 

q2‘(PD,pP)  (X)  =  q2(£;pD(X),PF(X))  (24) 

We  propose  to  use  the  golden-section  search  [21]  to  solve  (P4) .  This  requires 
that  q2*(PD,PF)  be  unimodal.  The  following  conditions  are  sufficient: 

1.  q2  is  strictly  convex 

2.  Pp  Pjj  and  Pp  ^  PF  ■  )  q2CS7Pp»PF)  q2C£?PjyPF) 

3.  PD(X)is  convex  n  and  PF{X)  is  convex  [J 

This  follows  since  for  every  X2,X2  e  A,  X2/X2, 

a  q2*(PD,PF)(X1)  +  (l-«)q2,(PD,PF!(X2) 

<  q2lS;  aPDtX^+d-ajPotXj),  «pf(x1)+(1-«)pf(X2)] 
iq2[s»  pD(«x1+(i-«)x2),  pF(«x1+(i-«)x2)] 

-  q2,(PD,PF)  («x1+(l-«)x2)  V  «  e  (0,1)  (25) 


>  *  — ;  p'lr  •*  - 

id Hat 
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and  so  q2*  (Pp»Pp)  *-s  strictlY  convex  A  and  thus  unimodal.  If  X=Pp,  then 
since  PD(X)  is  strictly  convex  A ,  the  conditions  simplify  to 

1.  q2  is  convex  A 

2.  For  every  Pp,  q2 (£;',Pp)  is  monotone  strictly  increasing 

and  similarly  for  x=PD.  These  conditions  seem  reasonable,  although  we  do  not 
give  proofs.  We  can  also  verify  that  q2* (PD,Pp)  is  unimodal  numerically. 

3.2  Posterior  Threshold  Optimization 

We  want  to  solve  the  following  functional  minimization  problem: 

(P5)  minimize  E{  i2} 

over  (PD,PF) :  Y*  — »  R 


Proceeding  as  in  the  prior  case  this  problem  is  reduced  to  the  pointwise 
minimization:'4. 


(P6)  for  any  given  Yk, 
minimize  tr  £^1^ 
over  (PD,PP)  6  R 


2Note  that  the  conditioning  here  is  actually  on  a  function  of  IPrjfPp) »  i.e., 
{<(PD,Pp),Yk(PD,PF)>:  (PD,Pp)  e  g},  as  opposed  to  a  value,  say  Y*(FD,Pp). 
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Let  (Pp,Pp)  be  the  optimal  value  of  (PD,Pp)  at  Yk  from  (P6) . 


Then 


(P&Pj)  (Yk)  =  (pJ,pJ) 


(26) 


is  optimal  for  (P5) .  Now 


m 

^k|k  =  tr  £k|k-l  -  t  E^jrtr(iikSk^)-|S£k2j|2] 

m 

+  E#i#j[aWti)'(Hi^j)]  }  (27) 

i/j»l 


Thus  tr  £^1^  will  be  minimized  when  the  quantity  in  {*}  is  maximized.  Hence 
our  problem  becomes: 


(P7)  for  any  given  Yk, 


maximize 


m 

Z#j(tr(HkSkH^)-|Hk2j|2] 


m 

+  Z^i^jKHk^i)'^)] 

i,j=l 

over  (PD,Pp)  e  R 


We  now  propose  an  algorithm  for  solving  (P7) .  First,  we  parameterize  PD,Pp  as 
in  the  prior  approach.  However,  instead  of  convexity  constraints,  we  require 
that  PD(X)  (or  equivalently,  Pp(*))  be  monotone  on  A  (for  convenience, 
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assume  Fj^x)  monotone  increasing  on  A  =  [*min'*maxl  )  •  'I^en  there  exist 
Xi,...,Xm*  such  that 

Xjnin  *  *0  <  Xi  <  ***  <  Xjjj*  <  -  Xj^  (28) 


and 


m=i  for  xe(Xi,Xi+i),  i=0,l,...,m* 


(29) 


(see  Figure  22) .  This  divides  (P7)  into  m*+l  local  maximizations  of  the  form 


m 


maximize  j  [tr  (Hj^H^)  -  iHfeStj  I  j 

j-1 

m 

+  S^jiawi)'^)] 
if  3*1 

over  x  6  [XifXi+1) 


which  we  may  write  as 


(P8)  maximize 


f — — 

Lb(X)+ 


)+a  (b(X)+a) 

over  x  e  [Xj.,Xi+1) 


■] 


I 

L 


where 
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b(X)  =  (2«)^2(VG/cMgM7c)PF(X)t(l-PD(>i)PG)/Pb(X)] 
m 

a  =  S  exp{-2'Sj^12j/2} 

j=l 

m 

dl  =  E  exp{-2’^1£j/2}  [ tr  (Hk^)  -  iHfcSj  1 2] 
m 

d2  =  E  exp{-£^12i/2}exp{£^1£j/2}  [  (H^)  *  (H^ )  ] 
i,j=l 


Note  that  if  and  d2  have  the  same  sign  (say  +) ,  then  (P8)  simplifies  to 


(P8')  minimize  Pp(X)  [(^(XJP^/PjjCX)] 
over  X  e  [*i,Xi+1) 


In  practice,  X1,...,Xm*  would  not  be  known  exactly  and  [X^,X^+^)  would  be 
approximated  from  within  by  a  closed  interval.  Finally,  we  note  that  this 
posterior  procedure  is  very  attractive  when  the  local  optimizations  [ (P8)  or 
(P8 ' ) ]  can  be  done  analytically. 

3.3  Simulation  Notes  and  Results 


The  simulation  of  the  prior  and  posterior  adaptive  threshold  optimization 
algorithms  involve  some  subtleties  which  we  shall  now  discuss.  First  note 
that  for  purposes  of  comparison  and  computational  efficiency  we  would  like  to 
use  the  same  Monte  Carlo  data  generated  for  the  fixed  (PD,Pp)  simulations  in 
Section  2.  In  the  prior  case,  for  every  Monte  Carlo  trial  there  will  in 
general  be  an  optimal  (PD,PF)  which  is  different  from  the  fixed  (PD,Pp) . 
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Thus,  for  every  Monte  Carlo  trial  we  need  to  evaluate  a  single  realization  of 
Yk  at  two  different  values  of  (PD,PF) .  This  can  be  done  in  two  ways: 

1.  Resetting  the  seed  in  the  generation  of  Yk 

k 

2.  Using  an  intermediate  result  from  a  single  generation  of  Y  (and 
saving  computation) 

Since  we  are  ultimately  concerned  with  Bernoulli  and  Poisson  randan  variates, 
approach  2  is  feasible  (we  emit  the  details) .  In  the  posterior  case,  we 
require  (and  not  just  for  purposes  of  comparison  and  computational  efficiency) 
a  single  realization  of  Y*  as  a  function  of  PD,PF,  or  sane  other  x.  It  turns 
out  that  approach  2  is  not  only  feasible,  but  exactly  specifies  the  values 
X-jy . . .  fXjn*  discussed  above.  Unfortunately,  time  was  not  available  to  run 
these  simulations. 


3.4  Conclusions 

We  have  carefully  posed  the  prior  and  posterior  detection  threshold 
optimization  problems,  and  have  given  algorithms  for  their  solution.  The 
prior  algorithm  performs  a  single  line  search  to  optimize  q2,  which  can  be 
evaluated  by  a  look-up  procedure.  The  posterior  algorithm  performs  multiple 
line  searches  which  can  sometimes  be  done  analytically.  Certain  subtleties  in 
the  simulation  of  these  algorithms  are  pointed  out,  and  the  simulations  turn 
out  to  be  neither  as  difficult  nor  as  expensive  as  originally  thought. 
Unfortunately,  time  did  not  permit  testing  on  either  real  or  simulated  data. 
It  is  expected  that  both  schemes  will  outperform  any  time- independent 
optimization.  It  is  also  expected  that  posterior  will  do  better  than  prior, 
although  it  may  be  more  expensive  computationally.  The  posterior  algorithm  is 
most  attractive  when  the  local  minimizations  can  be  done  analytically. 
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4.  GATE  SIZE  GPTIMIZATICN  AND  REVISE)  FDAF  DERIVATION 


Gate  size,  or  equivalently  Pq,  is  a  parameter  in  the  PDAF  equations.  Hence, 
Pq  could  be  subjected  to  the  same  analysis  given  PD'PF  above.  Specifically, 
we  could  obtain  a  measure  of  PDAF  performance  as  a  function  of  Pq  as  in 
Section  2,  and  then  optimize  over  Pq  using  similar  schemes  to  those  in  Section 
3.  However,  we  choose  not  to  pursue  these  directions  because: 

1.  PDAF  performance  is  approximately  constant  for  large  enough  gate 
size  (g  2  4  or  PG  2  .99) 

2.  PDAF  performance  is  approximately  monotonically  increasing  with  Pq, 
and  in  any  case,  joint  optimization  over  PD'PF  and  Pq  is 
considerably  more  complicated  than  just  over  PD,PF* 

For  these  reasons,  we  fix  Pq  in  a  time-independent  way  (g  =  4  or  Pq  -  1  for 
this  problem)  and  optimize  only  over  PD,Pp  as  in  Section  3.  However,  while 
considering  the  gate  size  issue,  we  observed  that  the  PDAF  exhibited  certain 
counter-intuitive  behavior.  We  next  address  this  problem. 


4.1  Gate  Size  in  PDAF  Derivation 


Consider  the  case  where  fl  =  I,  m  =  0,  PD  =  1,  ^^-1  »  £  and  Pq  ~  1  (g  2  4) . 
Since  m  =  0,  the  PDAF  sets  j^|k  =  £k|k-l  30(3  P^|k  =  £k|k-l*  But  w®  typically 
have 


16  <  Utklk-l4^3  'S^tfklk-r^  =  *k|k'*kjk-l*k|k  <30) 

which  implies  that  »  P^lk-l'  a  contradiction. 

Before  examining  the  PDAF  derivation  to  reconcile  this  contradiction,  we  give 
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a  simple  example  which  clearly  shews  hew  is  P°ssible  f°r  a 

Gaussian  problem  with  finite  gate.  Consider  the  problem  of  estimating  a 
scalar  x  -  N(0,p)  from  an  observation  y  given  by 


y  *  x+v 


(31) 


where  v  -  N(0,r)  and  x,v  are  independent.  It  is  well  known  that 


var(x|y)  ^  var(x) 


(32) 


However,  for  p  »  r  and  g  2  4,  we  have 


var (x|  lyl>gv^fr)  =  var(x|  |x|>  g^)  »  var(x) 


(33) 


The  point  is  that  conditioning  on  an  event  as  opposed  to  a  jointly  Gaussian 
random  variable  can  increase  the  conditional  variance.  It  should  be  clear 
that  this  simple  problem  has  all  the  relevant  features  of  the  case  described 
above. 

Turning  to  the  PDAF  derivation,  it  eventually  became  cle;  r  that  gate  size  was 
not  dealt  with  in  a  consistent  manner  or,  alternatively,  that  certain  unstated 
assumptions  were  made.  Me  next  show  how  to  update  the  conditional  mean  and 
covariance  when  gate  size  is  explicitly  accounted  for  and  also  give 
assumptions  for  which  the  PDAF  yields  the  same  results.  In  the  sequel,  we  use 
notation  as  developed  in  Section  1,  except  as  noted. 

Me  start  by  partitioning  the  event  Xg  into  the  events 
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X0£:  correct  measurement  is  not  detected 

X0gs  correct  measurement  is  detected  but  is  not  in  the  gate 


Note  that 


PtXgdiY1*}  =  PtXgdix^Y^pcxai^} 

=  PiXgdlXaJPtXgi^}  =  -^-PCXgiY*} 

1_PDPG 

and  similarly 

Ptx»IYk> 

To  compute 

*k|k  *  EUjJY15} 

*k|k  “  E**k|k*k|k|Yk} 


we  need  pU^lY*) .  But 


PlSklYk)  »  P{X0dlYk}p(ik|X0d,Yk)  +  PfXBgi^IpO^lXgg,^) 


m 


+  EPCXilY^p^IXirY*) 
i«l 


(34) 


(35) 


(36) 


(37) 


45 


Report  No.  5249 


Bolt  Beranek  and  Newman  Inc 


Previous  computation  of  P{Xj_|Yk},  i=0,l, . . . ,m,  was  correct,  so  using  standard 
notation  we  have 

pUklYk)  =  +  #0gP(*klX0g,Yk) 

m 

+  £  #iP(iSjclXi,Yk)  (38) 

i*l 

where 


*0d 


1-Pr 


i-Vo 


*0 


*0g  = 


*"Ws 


*0 


(39) 


Also, 


p(*k  1  ~  Ntfk|k-1'  ^k|k-l) 

P(*k  I  Xit  Yk)  -  fk|k-l«)'  i*l,...,m 


(40) 


Let  Gk  be  the  gate  at  the  k-th  iteration,  i.e. , 


Gjj  *  {£:  St'Sj^S  i  g2} 

! 


(41) 


From  Bayes  rule,  we  have 
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P(*k  I  *0q,Yk)  =  PUfc  I  Xflg,Yk'1,ni)  =  p(*k  I  X0g,Yk"1) 

■  P&k  1  ^  l) 

=  PtUSfcik-l+Sk^k  •  ^k^^lPtXk  •  ^k"1)/(1-pG) 

-  P^t  U2k  | k-l+%)  'Sk1  (EXk | k-l+iZk>  >  g 2}p(JSklYk"1)/(1-PG) 


where 


P&k)  -  N(0,B) 

pUklYk_1)  -  N(ik|k-!,  £k|k-l) 


Note  that  P„  { • }  is  in  general  very  hard  to  evaluate  (even  numerically) . 

*k 

At  this  point,  we  have  expressed  p(2Lk  lY*4)  in  terms  of  known  quantities.  We 
now  use  pC^lY*5)  ho  generate  Akjk  ^k|k*  P^rsh  consider  &k|k*  ^  have 


4ik 

”  ^0d^k|k-l  +  ^0gj*akP^iklx0g'Yk)^xk 


+  E^itiklk-l+W 

i*l 


JikPUkiXegr^xask  -  4ik-i 


(44) 
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or  equivalently 

J^lk-lPj^^lk-l^lk-l^Bg'^^ktk-l  *  4 

This  key  result  can  be  proved  as  follows.  It  is  sufficient  to  show 
every 

PaJt^k|k-l+^k|k-l^  =  (  ~^k  |  k-l+^k  |  k-1  ) 

and 

pXk { USk  |  k-i4^)  'Sk1  (B*k  |  k-14^  >  ^  * 

=  p^if-ifikik-i-^k)  '^(-flXkik-i-^k)  >g2> 

The  first  equation  follows  immediately  since 

pUk)  -  N^ik-i/Pklk-l^ 

and  the  second  follows  after  using 

P(Xk)  -  N(!,E) 


(45) 

that  for 

(46) 

(47) 

(48) 

(49) 
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(we  emit  the  details) .  Thus  we  have 


m 


4|k  =  *0d4|k-l  +  *0g4|k-l  +  2^*i(4|k-rfiik2i) 


=  4|k-l+4S  (50) 

which  is  the  samp  conditional  mean  update  as  the  PDAF.  Next  consider  £^1^. 
emitting  the  details,  we  have 


4|k  =j"4|k4|kP4lY 

=  4|k-l  +  ^0g(^|k-l-£k|k-l^  "  f1"*  0*444 
m 

+  Mk(E^j5tj23-22,)H4 

where 


(51) 


^lk-1  ®  ^4 1 k-l4  j k-lP 4 1 *0g' ^  3&k  (52> 

Comparing  this  recursion  with  the  stochastic  Riccati  equation,  we  see  there  is 
an  additional  term  *0g(  £fc|k-l-£k|k-l) '  increases  the  conditional 
covariance  update.  In  view  of  our  earlier  remarks,  this  is  not  unexpected. 

In  general,  cannot  be  expressed  analytically  in  terms  of  Pjtik-i.  It 
follows  that  if  no  approximation  is  made,  the  corresponding  filter  is 
infinite-dimensioned.  If  we  assume  that 
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PfeklXfl,**)  =  pkfclY^1)  (53) 

then  we  get  the  FDAF.  More  importantly,  the  FDAF  is  a  good  approximation  if 
and  only  if  *ggtr(  | k— l~^k | k— 1 )  is  suitably  small.  In  general,  0^ 

decreases  with  increasing  gate  size,  while  tr(pJ|k.1-jjl|k_1)  increases  with 
increasing  gate  size.  Consequently,  to  justify  using  the  FDAF,  gate  size  must 
be  chosen  large  enough  so  that  *0gtr  is  suitably  small,  and 

not  just  l-*0g  =  PG  =  1.  We  next  develop  a  method  to  compute 
approximately. 


4.2  Approximate  Computation  of 


We  start  by  approximately  computing  if  ^  for  two  extreme  cases.  We  have 

^c|k-l  =  E^k|k-A|k-1  I  *0g^} 

=  EUk|k-A|k-ll  ^Ik-l*^  '^^klk-l^^2'  Yk"1} 


r 

E^|k-l^k|k-l  I  ^k'B  1}  ■  £k|)c-i> 


BEkik-ifl'  «R 


(54) 


El*k\*-l*k\k-l  I  «k|k-lfl,(flEk|k-lfl')"1i2k|k-l>g2'  Y*"1}' 

® kik-ia' »» 


In  order  to  evaluate  the  »  £  case,  we  are  faced  with  the  problem  of 

computing  the  conditional  covariance  of  a  Gaussian  random  variable 
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given  that  jj^lk-l  lies  outside  a  hyperellipsoid  centered  at  its  conditional 
mean  (-£) .  Furthermore,  the  principal  axes  do  not  in  general  coincide  with 
the  eigenvectors  of  |k_^  (see  Figure  23a).  We  proceed  in  four  steps: 

1.  Perform  a  nonsingular  linear  transformation  on  which 

transforms  to  1  (see  Figure  23b). 

2.  Replace  the  hyperellipsoid  in  the  new  coordinate  system  with  a 
rectangular  approximation  (see  Figure  23c) . 

3.  Compute  the  conditional  covariance  of  2k|k-l  given  that  it  lies 
outside  the  rectangle. 

4.  Perform  the  inverse  transformation  of  step  1  to  get  an  approximation 
to  ^|k-l* 

We  omit  the  details  and  give  the  final  result: 


£g|k_l  =  Q' A1/2Q^ | fc-iQ' A172^  when  iffiklk-lH'  »  fi  (55) 


where 

Q  =  •••  %]  are  orthonormal  eigenvectors  of  |k— 1 

A  =  diag[X2  ...  Xn]  are  the  corresponding  eigenvalues 
1  *  A1/2Q'fl'  OSEk  | fc-ifi' ) _1a0A1/2 
3.  ■  [^2  ...  3n]  are  orthononnal  eigenvectors  of  Z 
A  =  diag[X2  ...  Xn]  are  the  corresponding  eigenvalues 

1-1  (g2^1  /2r3/2)  1-1  (g2^1  /2,3/2) 

■  9  •  •  •  9  ..  ■  ■  ■ 

.  2#(gRj)  2»(gK£) 


S|k-1  *  ^ 
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I(x,p)  =  U/rCp)] 


f 


VV1  du 


(incomplete  gamma  function) 


♦  (x)  S  [i/  2*] 


I 


“-u2/2 

x 


du 


(error  function) 


Since  I  is  in  general  not  invertible,  5  is  formed  by  augmenting  the 
orthonormal  eigenvectors  of  Z  (A  and  (k— 1  are  chan9ed  appropriately).  We 
note  that  IMSL  routines  exist  for  the  computation  of  I(x,p)  and  0  (x) . 

In  the  general  case  where  |  ^_iii '  r  E,  ^|k-l  is  approximated  by 

interpolating  between  the  two  extreme  cases  in  a  reasonable  manner. 


4.3  Conclusions 


We  started  by  pointing  out  that  gate  size  or,  equivalently  PG,  could  be 
subjected  to  the  same  analysis  given  PD,Pp  in  Sections  2  and  3.  Specifically, 
we  could  obtain  a  measure  of  performance  as  a  function  of  PQ  and  then  optimize 
over  Pq.  For  various  reasons  it  appears  uninteresting  and/or  unproductive  to 
continue  in  this  direction.  However,  in  considering  the  gate  size  issue,  it 
was  observed  that  the  PDAF  exhibited  certain  counter-intuitive  behavior.  Upon 
examination  of  the  PDAF  derivation  it  became  clear  that  gate  size  was  not 
dealt  with  in  a  consistent  manner,  or  alternatively,  that  certain  unstated 
assumptions  were  made.  We  have  shown  how  to  update  the  conditional  mean  and 
covariance  when  gate  size  is  explicitly  accounted  for  and  have  also  given 
assumptions  under  which  the  PDAF  yields  the  same  results,  it  turns  out  that 
the  conditional  mean  update  is  the  same  as  in  the  PDAF,  but  the  conditional 
covariance  update  is  increased  by  an  additional  term.  A  method  is  developed 
to  approximate  this  additional  term  by  considering  the  extreme  cases  where 
process  noise  »  measurement  noise  and  conversely.  The  derivations  in  this 
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section  are  important  for  theoretical  completeness,  and  may  also  be  of 
practical  importance  under  certain  operating  conditions  (e.g. ,  if  one  selects 
a  moderate  gate  size  for  computational  reasons).  No  simulations  were  run  for 
lack  of  time  and  we  leave  this  as  a  future  task. 
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